clc
clear
close all
SampFreq = 10000;
t = 0:1/SampFreq:1;
IFrot = 20*t.^2-40*t+45; % Instantaneous shaft rotating frequency
phirot = 2*pi*(20/3*t.^3 - 20*t.^2 + 45*t); %Phase of the signal
FCC = 5;% Fault characteristic co-efficient
IFFCF = FCC * IFrot; %Instantaneous fault characteristic frequency
phiFCF = FCC * phirot; %Phase of fault characteristic frequency
FR = 2000; %Resonance frequency
alpha = 0.9; 
lamda = 0; 
beta = 1200; % Damping co-efficient
amp = 10 + alpha*cos(phirot);
[~, index] = impultime(phiFCF,t);
comimpulse = exp(-beta*t).*sin(2*pi*FR*t);% Damping of the impulses
impulse = zeros(1,length(t));
for ii = 1:length(index)
    subimpulse = zeros(1,length(t));
    subimpulse(index(ii):end) = comimpulse(1:end-index(ii)+1);
    impulse = impulse + (1 + lamda*IFrot(index(ii))).*subimpulse;
end
impulse2 = amp.*impulse;

plot(t,impulse2);
xlabel('Time/s')
ylabel('Amplitude')
title('Damping co-efficient is 3600')